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ANISOTROPY  IN  REACTING  COMPRESSIBLE 
TURBULENT  SHEAR  FLOW 


D.  LIVESCU  AND  C.K.  MADNIA 

Department  of  Mechanical  and  Aerospace  Engineering, 

State  University  of  New  York  at  Buffalo, 

Buffalo,  NY  14260 


Abstract. 

The  structure  of  reacting  compressible  homogeneous  turbulent  shear  flow  is 
examined  via  data  generated  by  direct  numerical  simulations  (DNS).  The  reaction 
is  modeled  as  one-step,  exothermic,  irreversible,  and  Arrhenius  type.  The  effects 
of  reaction  are  studied  by  decomposing  the  velocity  field  into  its  solenoidal  and  di- 
latational  parts.  The  heat  of  reaction  significantly  enhances  the  dilatational  kinetic 
energy,  while  it  reduces  its  solenoidal  counterpart.  The  solenoidal  and  dilatational 
large  scales  (kinetic  energy)  and  small  scales  (viscous  dissipation)  anisotropies  are 
investigated  for  cases  with  different  initial  turbulent  Mach  numbers  and  Reynolds 
numbers. 


1.  Introduction 

Homogeneous  turbulent  shear  flow  represents  one  of  the  simplest  anisotropic  flows 
and  its  study  can  reveal  important  aspects  of  the  structure  and  production  of  the 
turbulent  fluctuations.  Moreover,  the  recent  high  Reynolds  number  experiment  of 
incompressible  turbulent  shear  flow  of  Shen  and  Warhaft  (2000)  indicates  that  the 
higher  order  moments  of  the  velocity  field  are  not  consistent  with  the  postulate 
of  local  isotropy.  Although  DNS  at  high  Reynolds  numbers  are  not  feasible  yet, 
persistent  anisotropy  of  the  higher  order  moments  for  incompressible  shear  flow 
at  moderate  Reynolds  numbers  has  been  observed  earlier  in  the  numerical  results 
of  Puinir  (1996).  Compared  to  incompressible  turbulent  shear  flow,  the  compress- 
ible turbulent  shear  flow  is  much  less  investigated  computationally,  although  the 
recent  studies  of  Blaisdell,  Mansour  & Reynolds  (1993),  Sarkar  (1995)  and  Hamba 
(1999)  have  improved  our  understanding  of  the  influence  of  compressibility  on  the 
structure  and  development  of  the  flow.  The  heat  release  influence  on  the  turbu- 
lent shear  flow  development  has  been  studied  by  Livescu,  Jaberi  & Madnia  (2001) 
(hereinafter,  referred  to  as  LJM01).  However,  the  reacting  flows  usually  involve  a 
large  amount  of  heat  release  and  thus  are  intimately  related  to  the  compressibility 
of  the  flow.  Furthermore,  the  combined  effects  of  heat  release  and  compressibility 
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on  the  structure  of  the  turbulent  shear  flow  have  not  been  studied.  The  main  objec- 
tive of  this  study  is  to  examine  the  influence  of  heat  release  and  compressibility  on 
the  large  scale  and  small  scale  anisotropies  in  reacting  compressible  homogeneous 
turbulent  shear  flow. 


2.  Numerical  Methodology  and  DNS  Parameters 

In  order  to  assess  the  coupled  influence  of  compressibility  and  heat  release  on  the 
structure  of  the  reacting  compressible  shear  flow,  DNS  of  homogeneous  sheared 
turbulence  are  performed  under  reacting  (heat-releasing)  and  nonreacting  con- 
ditions. The  continuity,  momentum,  energy  and  species  mass  fractions  transport 
equations  are  solved  using  the  spectral  collocation  method.  The  coordinate  system 
is  moving  with  the  mean  velocity  which  makes  necessary  a periodic  remeshing  of 
the  grid  in  order  to  avoid  errors  associated  with  highly  skewed  grids.  The  remesh- 
ing procedure  is  carried  out  in  wavenumber  space  with  explicit  removal  of  the 
modes  affected  by  aliasing.  The  aliasing  errors  associated  with  the  evaluation  of 
the  nonlinear  terms  are  controlled  by  using  the  skew-symmetric  form  of  the  con- 
vective terms.  In  the  reacting  cases,  the  chemical  reaction  is  modeled  by  a single 
step,  exothermic,  irreversible,  Arrhenius-type  reaction.  The  viscosity  varies  with 
the  temperature  according  with  a power  law  and  Le  = 1 in  all  cases. 

The  velocity  fluctuations  are  initialized  as  a random  solenoidal,  3-D  field  with 
Gaussian  spectral  density  function  (with  the  peak  at  kov  = 10)  and  unity  rms. 
The  initial  pressure  fluctuations  are  evaluated  from  a Poisson  equation  (except  for 
case  3 for  which  they  are  set  to  zero)  and  the  initial  density  field  has  the  mean 
value  equal  to  1 and  no  fluctuations.  For  the  reacting  cases  the  scalar  fields  are 
initialized  as  “random  blobs”,  with  double-delta  PDFs  (Overholt  & Pope,  1996). 
Details  can  be  found  in  LJM01.  The  variables  are  time  advanced  in  physical  space 
using  a second  order  accurate  Adams-Bashforth  scheme.  The  code  developed  for 
this  work  is  based  on  a fully  parallel  algorithm  and  uses  the  standard  Message 
Passing  Interface  (MPI). 

The  cases  considered  in  this  study,  summarized  in  table  1,  have  different  values 
for  the  initial  turbulent  Mach  number,  Mto , and  initial  Reynolds  number,  Re\0 . All 
cases  have  the  same  value  of  the  mean  shear  rate  S = du\  /dx2  = 5.1  (Hi  is  the  mass 
average  of  the  velocity),  which  is  in  the  range  dominated  by  nonlinear  effects.  For 
each  case  two  simulations,  nonreacting  and  reacting,  are  performed.  All  reacting 
cases  have  the  values  of  the  heat  release  parameter  Ce  — 1.44,  computational 
Damkohler  number  Da  = 1100,  and  Zeldovich  number  Ze  = 8.  For  an  initial 
mean  temperature  in  the  range  of  temperatures  prior  to  ignition  in  an  IC  engine, 
the  values  for  the  pre-exponential  constant  of  reaction  Kf,  the  activation  energy 
Ea,  and  the  heat  of  reaction  —H°,  are  calculated  and  they  are  in  the  range  of 
values  for  the  elementary  reactions  for  hydrocarbon  combustion. 

The  range  of  Mach  numbers  considered  extends  from  the  nearly  incompressible 
case  1 to  the  high  Mach  number  case  3 which  is  at  the  upper  limit  for  the  nu- 
merical method  considered.  The  range  of  Re\0  extends  to  values  higher  than  any 
considered  in  previous  studies  of  compressible  homogeneous  turbulent  shear  flow. 
The  simulations  were  monitored  to  ensure  that  the  integral  scales  remain  much 
smaller  than  the  box  size,  all  scales  of  motion  are  fully  resolved  and  no  significant 
shocklet  regions  develop.  Thus,  all  simulations  were  stopped  between  St  = 18  and 
20,  except  the  nonreacting  and  reacting  case  5 which  were  stopped  at  St  = 10.5 
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Case  # 

Mto 

Re  x0 

Grid  size 

1 

0.1 

21 

1283 

2 

0.3 

21 

1283 

3 

0.6 

21 

1283 

4 

0.3 

14 

643 

5 

0.3 

50 

2563 

2+ 

0.3 

21 

1283 

5+ 

0.3 

50 

2563 

TABLE  1.  Parameters  for  the  DNS  cases.  + In 
these  cases  5 is  set  to  zero  at  St  = 4. 


and  St  = 12.5,  respectively,  before  strong  shocklets  can  be  formed  in  the  flow. 
Cases  with  different  values  of  the  reaction  parameters  and  mean  shear  rate  have 
also  been  examined.  Although  the  results  for  these  cases  are  not  shown  here,  they 
agree  qualitatively  with  the  results  presented  in  this  paper. 


3.  Results 

Due  to  the  presence  of  the  mean  velocity  gradient,  the  velocity  field,  initially 
isotropic,  develops  to  become  anisotropic  with  K\  > K3  > K3,  where  K\  is  the 
kinetic  energy  in  the  direction  of  the  mean  velocity  (the  only  component  of  the 
kinetic  energy  with  a production  term  in  its  transport  equation)  and  K3  is  the 
kinetic  energy  in  the  direction  of  the  shear.  In  order  to  examine  quantitatively  the 
degree  of  anisotropy  among  the  kinetic  energy  components,  the  Reynolds  stress 
anisotropy  tensor  is  considered: 


< pu'i'u'j  > 1 

< P«  > 3 ij 


(1) 


In  agreement  with  the  previous  nonreacting  simulations  (Sarkar  1995;  Hamba 
1999)  our  results  show  that  the  anisotropy  among  the  normal  stresses  increases 
with  increasing  Mto.  However,  as  figure  1 shows,  for  the  reacting  cases  the  normal 
stresses  anisotropy  decreases  during  the  time  when  the  reaction  is  important  (3  < 
St  < 8).  Similar  results  are  also  obtained  for  reacting  mixing  layers  (Luo  1999). 
For  the  case  of  reacting  homogeneous  shear  flow,  it  is  shown  in  LJM01  that  this 
decrease  can  be  associated  with  two  effects,  a decrease  in  the  solenoidal  kinetic 
energy  in  the  direction  of  the  mean  velocity  due  an  enhanced  viscous  dissipation 
and  an  increase  in  the  dilatational  part  of  the  kinetic  energy  in  the  direction  of 
the  mean  shear.  However,  after  the  reaction  is  completed  (St  > 10)  the  normal 
stresses  anisotropy  becomes  larger  than  in  the  nonreacting  case.  This  behavior  is 
related  to  reduced  levels  of  the  pressure  strain  terms,  responsible  for  the  energy 
redistribution  among  the  kinetic  energy  components.  As  the  value  of  the  shear 
rate  5 is  set  to  zero,  the  flow  starts  to  return  to  isotropy  (figure  1),  although  some 
small  persistent  anisotropy  can  be  observed  for  the  time  range  examined. 
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As  explained  by  LJM01,  the  solenoidal  and  dilatational  parts  of  the  kinetic 
energy  are  influenced  differently  by  the  reaction.  Therefore,  in  order  to  understand 
the  heat  release  effects,  it  is  more  useful  to  study  separately  the  anisotropy  among 
the  dilatational  and  solenoidal  components  of  the  normal  Reynolds  stresses.  We 
define  a measure  of  the  anisotropy  of  the  normal  stresses  as: 

«=[(^i)2  + (^2)2  + (%)2]1/2  (2) 

where  bfj  =<  pu"au'ja  > / < > —1/3%,  and  a = s or  d represents  the 

solenoidal  or  dilatational  part,  respectively.  Helmholtz  decomposition  is  used  to 
split  the  normal  Reynolds  stresses  into  their  solenoidal  and  dilatational  compo- 
nents. In  the  extreme  case  where  6fx  = 2/3  and  622  = ^33  = -1/3,  = ^/2/3 

and  thus  0 < a%  < ^2/3. 

Figure  2 shows  that  aK  increases  as  Mt0  increases,  and  this  can  be  mainly 
associated  with  a decrease  in  the  solenoidal  pressure  strain  terms.  For  the  reacting 
cases,  as  shown  by  LJM01,  the  heat  release  influences  the  solenoidal  parts  of  the 
normal  Reynolds  stresses  mainly  through  the  viscous  terms.  Furthermore,  the  in- 
crease in  the  solenoidal  viscous  dissipation  components  is  such  that  the  anisotropy 
among  the  solenoidal  parts  of  the  normal  stresses  is  little  affected.  As  a result,  as 
figure  2 shows,  the  values  of  a*K  for  the  reacting  and  nonreacting  cases  are  close 
during  the  time  when  the  reaction  is  important.  At  later  times  asK  becomes  larger 
for  the  reacting  cases  compared  to  the  nonreacting  cases,  although  the  difference 
between  the  values  obtained  for  the  reacting  and  nonreacting  cases  decreases  with 
increasing  Mach  number. 

The  influence  of  the  initial  value  of  the  Reynolds  number  on  the  evolution  of  asK 
is  presented  in  figure  3(a).  As  expected,  a$K  decreases  with  increasing  Re\0  for  both 
nonreacting  and  reacting  cases.  After  S is  set  to  zero,  asK  decreases  initially,  then 
reaches  a quasi-steady  state.  This  persistent  (although  small)  anisotropy  decreases 
with  increasing  Re\0  and  it  is  discussed  in  more  detail  below. 

It  is  found  in  LJM01  that  the  explicit  dilatational  effects  occur  predominantly 
in  the  direction  of  the  mean  velocity  gradient  in  both  nonreacting  and  reacting 
cases.  The  role  of  the  reaction  is  to  significantly  amplify  these  effects,  during 
the  time  when  the  reaction  is  important.  Consequently,  as  figure  3(b)  shows,  the 
anisotropy  among  the  dilatational  components  of  the  normal  Reynolds  stresses 
increases  by  the  reaction.  By  comparing  the  results  presented  in  figures  3(a)  and 
(b)  it  can  be  seen  that  the  anisotropy  levels  remain  higher  for  a ^ than  for  ai- 
after  5 is  set  to  zero.  For  clarity,  the  results  obtained  for  the  reacting  cases  2^ 
and  5+  are  not  shown  in  figure  3(b)  but  they  are  consistent  with  the  discussion 
above.  The  Mt0  influence  on  the  evolution  of  adK  has  also  been  examined  but  no 
clear  trend  has  been  observed,  although  it  seems  that  the  behavior  at  small  and 
large  values  of  Mf0  is  different. 

The  quantity  a%  is  defined  based  on  the  kinetic  energy  components  and  thus 
it  can  be  associated  with  the  anisotropy  at  large  scales.  The  small  scale  anisotropy 
has  been  studied  in  detail  for  incompressible  flows  for  a large  range  of  Re \ and  it 
is  generally  agreed  that  the  second  order  moments  are  consistent  with  the  require- 
ments of  local  isotropy  postulate  (Saddoughi  & Veeravalli  1994;  Shen  & Warhaft 
2000;  Ferchichi  & Tavoularis  2000).  Much  less  is  known,  however,  about  the  small 
scale  anisotropy  in  compressible  flows  or  the  differences  between  the  small  scale 
behavior  of  the  solenoidal  and  dilatational  velocity  fields.  Although  the  values  of 
the  Reynolds  numbers  in  the  present  study  are  not  high  enough  to  address  the 
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issues  of  the  local  isotropy  postulate,  our  results  can  be  useful  for  modeling  pur- 
poses and  to  shed  some  light  on  the  physics  of  the  solenoidal  and  dilatational  fields 
in  turbulent  reacting  flows. 

Consistent  with  the  definitions  for  the  viscous  dissipation  terms  in  the  trans- 
port equations  for  the  solenoidal  and  dilatational  parts  of  the  normal  Reynolds 
stresses  provided  in  LJM01,  the  measures  of  the  anisotropies  of  the  solenoidal  and 
dilatational  dissipations  are  defined  as: 

<=[K!)2  + (42)2  + (^3)2]1/2  (3) 

where  dfj  = - 1 / , and  efj  and  efj  are  the  components  of  the  solenoidal 

and  dilatational  viscous  dissipation  tensors,  respectively.  Figure  4 shows  that  the 
solenoidal  dissipation  anisotropy  has  a similar  dependence  on  M(o  as  a Jo  for  non- 
reacting as  well  as  for  reacting  cases.  The  difference  between  the  late  time  values 
of  asK  and  af  for  the  reacting  and  nonreacting  cases  decreases  at  higher  values  of 
M,0 . In  order  to  examine  which  scales  of  motion  are  responsible  for  this  behavior, 
the  three-dimensional  power  spectra  of  the  solenoidal  velocity  are  examined  in 
figure  5 at  St  = 18.  It  can  be  seen  that  the  energy  levels  are  lower  at  all  scales 
at  higher  Mto  for  both  nonreacting  and  reacting  cases.  However,  for  the  reacting 
cases  the  small  and  intermediate  scales  are  mostly  affected.  Since  the  small  scales 
are  more  isotropic,  the  energy  decrease  at  small  scales  for  the  reacting  cases  can 
be  associated  with  an  increase  in  the  anisotropy  of  the  flow.  Moreover,  for  the 
reacting  cases  the  small  scales  energy  levels  decrease  more  at  lower  Mt0.  As  the 
Reynolds  number  increases,  a*  decreases  (figure  6a).  However,  it  remains  larger  for 
the  reacting  cases  than  for  the  nonreacting  cases.  For  the  highest  Re\0  examined, 
as(  becomes  very  small  after  S is  set  to  zero. 

The  dilatational  dissipation  anisotropy  becomes  larger  than  its  solenoidal  coun- 
terpart after  the  initial  development  time  for  both  reacting  and  nonreacting  cases 
(figure  6b).  Unlike  a®,  af  is  much  larger  for  the  reacting  cases  than  for  the  corre- 
sponding nonreacting  cases  during  the  time  when  the  reaction  is  important,  which 
is  consistent  with  the  discussion  above  about  the  amplification  of  the  dilatational 
effects  by  the  reaction  in  the  direction  of  the  mean  velocity  gradient.  After  the 
mean  shear  rate  is  set  to  zero,  even  for  the  case  with  the  highest  Re\0  examined, 
there  are  still  some  persistent  anisotropy  levels.  For  clarity,  the  results  obtained  for 
the  reacting  cases  2+  and  5+  are  not  shown  in  figure  6(b)  but  they  are  consistent 
with  the  discussion  above. 

For  reacting  decaying  isotropic  turbulence,  Jaberi,  Livescu  & Madnia  (2000) 
show  that  the  heat  release  increases  the  dilatational  kinetic  energy,  unlike  the 
solenoidal  energy,  at  all  scales.  Figure  7 shows  that  at  St  = 6,  which  is  close  to 
the  time  when  the  mean  reaction  rate  peaks  for  case  5,  the  dilatational  energy  is 
enhanced  at  all  scales  by  the  reaction.  As  the  reaction  is  completing  (St  > 10),  the 
contributions  from  the  production  and  pressure  dilatation  terms  to  the  increase 
in  the  dilatational  kinetic  energy  decrease  (LJM01).  Moreover,  due  to  the  higher 
temperature,  the  viscous  dissipative  effects  are  still  significant  and  the  energy  in 
the  small  dilatational  scales  becomes  less  than  in  the  nonreacting  case. 

In  order  to  understand  the  persistent  anisotropy  levels  observed  in  the  normal 
Reynolds  stresses  after  setting  S to  zero,  the  three-dimensional  power  spectra  of 
the  solenoidal  and  dilatational  velocities  for  the  nonreacting  case  5+  at  different 
times  are  examined  in  figure  8.  It  can  be  seen  that  the  solenoidal  kinetic  energy 
is  affected  mostly  at  intermediate  scales  after  the  shear  rate  is  set  to  zero.  Thus 
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the  largest  eddies,  which  are  more  anisotropic,  are  little  affected  during  the  initial 
transient  time  when  the  flow  adapts  to  the  new  conditions.  However,  the  inter- 
mediate dilatational  scales  do  not  decrease  their  energy  during  initial  transient 
time  after  setting  5 = 0,  and  the  largest  scales  even  increase  their  energy  levels. 
Therefore,  it  is  expected  that  some  persistent  dilatational  anisotropy  levels  should 
remain  in  the  flow  after  setting  5 = 0.  Although  these  levels  are  likely  to  vanish 
as  the  Reynolds  number  increases,  it  is  important  to  note  that  they  are  larger 
than  for  the  solenoidal  motions,  suggesting  that  the  return  to  isotropy  could  be 
different  for  the  solenoidal  and  dilatational  fields. 


4.  Concluding  remarks 

DNS  are  conducted  of  compressible  homogeneous  turbulent  shear  flow  under  react- 
ing (heat  releasing)  and  nonreacting  conditions  to  examine  the  coupled  influence 
of  compressibility  and  heat  release  on  the  structure  of  the  flow.  For  the  reacting 
cases  the  chemical  reaction  is  modeled  as  one  step,  irreversible,  and  Arrhenius 
type.  The  parameters  considered  for  this  study  are  the  initial  values  of  the  turbu- 
lent Mach  number  and  the  Reynolds  number.  In  order  to  examine  the  return  to 
isotropy,  after  significant  anisotropy  levels  are  developed  in  the  flow,  the  value  of 
the  mean  shear  rate  is  set  to  zero  for  cases  with  different  values  of  Re\0 . 

In  order  to  study  the  changes  in  the  structure  of  the  flow  due  to  the  heat  of 
reaction  and  compressibility,  the  velocity  field  is  decomposed  into  solenoidal  and 
dilatational  parts  using  the  Helmholtz  decomposition.  Anisotropy  tensors  are  de- 
fined for  each  part  of  the  normal  Reynolds  stresses  and  also  for  the  solenoidal  and 
dilatational  parts  of  the  viscous  dissipation.  The  results  show  that  both  the  large 
and  small  scale  anisotropies  are  affected  by  compressibility  and  heat  release.  How- 
ever, the  small  solenoidal  scales  tend  to  be  more  isotropic  than  their  dilatational 
counterparts  for  the  range  of  Reynolds  numbers  examined.  This  tendency  is  also 
maintained  after  setting  the  shear  rate  to  zero. 
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Figure  1:  The  diagonal  components  o!  the  anisotropy 
tensor  lor  case  2 (thin  lines)  and  2*  (thick  lines). 


Figure  2:  Mach  number  effect  on  the  evolution  of  aK’  for 
the  nonreacting  (thin  lines)  and  reacting  (thick  lines)  cases. 
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Figure  3:  Re  influence  on  the  evolution  of  a)  aK‘  and  b)  aK  for  the  nonreacting  (thin  lines)  and  reacting  (thick  lines) 
cases. 


